Tempo and Mode of Evolution in the Tangled Nature Model 
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We study the Tangled Nature model of macro evolution and demonstrate that the co-evolutionary 
dynamics produces an increasingly correlated core of well occupied types. At the same time the 
entire configuration of types becomes increasing de-correlated. This finding is related to ecosystem 
evolution. The systems level dynamics of the model is subordinated to intermittent transitions 
between meta-stable states. We improve on previous studies of the statistics of the transition times 
and show that the fluctuations in the offspring probability decreases with number of transitions. 
The longtime adaptation, as seen by an increasing population size is demonstrated to be related to 
the convexity of the offspring probability. We explain how the models behaviour is a mathematical 
reflection of Darwin's concept of adaptation of profitable variations. 
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PACS numbers: 

I. MOTIVATION 

In his 1944 book [l|, Simpson's expressed aim was 
to reconcile the fruits of paleontological research with 
the Darwinian synthesis. The Darwinian synthesis, such 
as it was, consisted of Darwin's original ideas concern- 
ing competition for scarce resources and natural selec- 
tion married to more recent discoveries revealing the first 
genetic, then fundamentally molecular, mechanisms by 
which Darwin's ideas were realised. The question that 
Simpson wanted to answer was: is it possible to infer 
anything about the mode of evolution simply by observ- 
ing its tempo. 

The point, that has been made in different ways since 
Simpson's seminal work, is that, as Simpson has it, and 
as Gould has subsequently underlined |2j , what holds for 
ten rats over a hundred years may not be true, or may 
not be simply extrapolated, to a million rats over a billion 
years. In this paper, we take this point further, by ob- 
serving the tempo of evolution of a model of interacting 
populations of different species. We find that there are 
two separate dynamical regimes on which evolution takes 
place. In one of those, the tempo gradually decelerates 
due to the frustration in the interactions between differ- 
ent species. In the other, observable quantities fluctuate 
around stationary values. 

A related issue is whether evolutionary dynamics on 
all spatial scales, ranging from bacterial colonies con- 
fined in a flask (see eg. 3) to biological species can be 
understood, at least in an approximate fashion, using a 
small set of unifying theoretical principles. These prin- 
ciples would explain the emergent properties of evolu- 
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tionary systems, for example the decrease of the rate of 
evolutionary events on coarse-grained time scales in 
spite of the microscopic rate of change remaining con- 
stant. Clearly, evolution at the collective level is a non- 
stationary process, characterised by various trends. How 
are these trends related to the dynamics at the level of 
individuals? . 

An approach which has gained momentum due to 
the increasing speed and availability of computers, is to 
search for unifying principles by simulating simple evo- 
lutionary models with Monte Carlo techniques. These 
models blend a random element, the noise, with deter- 
ministic rules for the average rates of reproduction, mu- 
tation and death in a group or groups of individuals. 
Precursors to our approach include Kauffman's NK(C) 
model, which has been used primarily as a tool to study 
the way fitness landscapes affect the course of evolution 
Q . Our approach, while partially derived from the NKC 
model, is different in two ways: we have frequency de- 
pendent interaction, and we do not explicitly build in a 
optimisation. The fact that our model nevertheless ap- 
pears to optimise certain macroscopic parameters is a 
matter of some interest, and is discussed further below. 

There are many elements of the puzzle that we will 
simply neglect. Indeed the model we shall study does 
not even have the rudimentary innovation of sexual re- 
production. This approach, we should make clear, is lim- 
ited and cannot provide a complete answer to the ques- 
tion of how the microscopic and macroscopic are related. 
However, we show that even such a simple model as the 
one discussed can give rise to macroscopic evolution of 
a very different tempo from the individual dynamics on 
which it is based, and we show how this emergent macro- 
scopic behaviour comes about from the interaction of 
many individuals following Darwinian principles. It is 
the glassy nature of evolutionary dynamics - the exis- 
tence of quenched disorder giving rise to slow relaxation 
- that in our model produces the macroscopic tempo of 
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a different nature from that of the underlying dynamics. 



II. THE MICROSCOPIC LEVEL 

The Tangled Nature model of co-evolution [6j has al- 
ready been studied in several contexts. Its simplic- 
ity along with the rich complexity of its resulting be- 
haviour makes it a paradigmatic model for testing co- 
evolutionary ideas. The model retains the binary string 
genotype geometry found in previous approaches, but 
replaces their 'ad hoc' static fitness landscapes with a 
set of population dependent interactions between extant 
species, similar to the 'tangled' interactions of an eco- 
system. This simple addition creates a fascinating com- 
plex behaviour: the model retains the intermittent dy- 
namics found for adaptive walks in fitness landscapes [? 
], where periods of quasi-stasis which are long with re- 
spect to the timescale of the underlying, reproductive 
dynamics, are interspersed by transition events, termed 
quakes, which occur at a decreasing rate. The overall 
tempo of evolution is decelerating in the model, indicat- 
ing the possibility of an underlying process optimising 
a hidden variable. Importantly, the model introduces 
new, and biologically relevant, variables in the dynamics 
at the macroscopic level. As the system evolves from a 
'random' initial state in the presence of scarce external 
resources, the network of extant and interacting popula- 
tion changes, slowly, but radically over time, enabling the 
system to support an ever growing number of individuals 
in the face of a constant supply of external resources. 

Despite its simplicity, the model is able to reproduce 
the long time decrease reported in the overall macro- 
scopic extinction rate, the observed intermittent nature 
of macro-evolution, denoted punctuated equilibrium by 
Gould and Eldrdege [8|, the log- normal shape often ob- 
served for the species abundance distributions @, the 
power law relation often seen between area and the num- 
ber of different species number 10], the frame work of 
the model is also able to reproduce often reported expo- 
nential degree distributions of the network of species as 
well as the decreasing connectance with increasing species 
diversity that has attracted much observational and the- 
oretical interest [i~lj |. 

Individuals are defined by a binary string of length 
L, and are thus points in an L dimensional hypercube, 
whose elements, or nodes, will be denoted by Latin letters 
i and j. An a priori given fixed and random matrix Jy , 
specifies the strength and type of potential interaction 
between populations at nodes i and j. The evolutionary 
dynamics of the model, to be explained in detail below, 
leads to the occupancy of types evolving with time. Mu- 
tation allows new nodes to become inhabited and existing 
occupied types to go extinct. The configuration of occu- 
pied types changes accordingly with time. Although the 
dynamics at the level of individuals is smooth and al- 
ways happens at the same rate, at the systems level, i.e. 
the level of occupied configurations in type space, the 



dynamics is highly intermittent, and, as mentioned, de- 
celerating. The dynamics of the model is characterised 
by a succession of quasi stable states (or QESS: quasi 
evolutionary stable strategy). Each QESS is described 
by a network of sites which are populated and interact 
with each other. As the system undergoes the intermit- 
tent burst of changes to the occupancy one network of 
occupied sites is replaced by another. 

The long time evolution of the occupancy in type 
space, and all the properties that can be derived form the 
composition and interactions of the extant population, 
are subordinated to the transitions, or the quakes, as 
they are customarily denoted in record dynamics. This is 
because no essential change occurs during the metastable 
state. The occupancy of individual sites may fluctuate 
but, since essentially the same set of types are present 
during the metastable periods, no net drift can manifest 
itself. The most conspicuous overall effect of the suc- 
cession of QESS is a logarithmically slow increase in the 
overall population size. 

The evolution of the configuration in type space is non- 
stationary nature and the population grows logarithmi- 
cally with time as the configurations, through the adap- 
tive search that takes place during the quakes events, 
manage to become collectively better adapted. The 
adaptation consists in selecting networks of occupied sites 
in type space connected by increasingly mutualistic or 
beneficial interactions, which enables the temporal av- 
eraged net-reproduction rate to slightly overcome the 
killing rate. 

A time-step in the model is carried out as follows: an 
individual is picked at random. Its offspring probability 
is calculated by first calculating the value of 

j 

where N is the total population, rij the site popula- 
tions, C a strength parameter and /i the carrying capac- 
ity determining the resources available to the population. 
This value is then mapped on to the interval 0,1 using 
the function 



° ff l + exp(-tf)' 

which is the probability of reproduction in that time 
step. Reproduction is asexual: the parent is replaced 
by two copies of itself, with each offspring susceptible to 
mutations with probability p mu t per genome site. In the 
next step a randomly chosen individual is killed with con- 
stant probability pkui ■ Note that while a point mutation 
creates a new individual on a neighboring site of the hy- 
percube on which populatons reside, the interaction ma- 
trix generally involves populations located at arbitrary 
positons on the hypercube. 

Typically the population is initially located on a ran- 
domly chosen, single site. To start with, the population 
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decreases. Then, as reproduction becomes more likely, 
mutants begin to populate neighbouring sites. Soon 
a stable "ecosystem" is formed, a subset of genotypes 
whose mutual interactions form a (meta) stable dynami- 
cal system. This subsystem can be characterised by some 
key macroscopic properties - the average total number of 
individuals in the system, the average number of sites 
occupied, and the properties of the realised interaction 
network. 



III. THE MACROSCOPIC LEVEL 

A. Centre of mass 

Since the system is made up of weighted nodes (each 
with a vector position in genotype space) interacting with 
each other, we can define a quantity analogous to the 
centre of mass in classical mechanics. We define this as 
the species occupancy weighted average position of the 
system in genotype space, A, where 

A= jvE rijSi W 

j 

where Sj denotes the position of vertex number j = 
1,...,2 L in the L dimensional hypercube of genotypes, 
fij is the occupancy of species j located at vertex Sj and 
N is the total population: (N — J^j n j)- 

During meta-stable periods this quantity remains prac- 
tically constant, since it is only slightly modified by pop- 
ulation fluctuations, which chiefly occur on the periph- 
ery of the system due to mutations. In contrast, when a 
quake occurs the centre of mass undergoes a major dis- 
placement. While it may be theoretically possible that 
a quake would occur which does not shift the centre of 
mass perceptibly, this is highly unlikely since any nearby 
stable configuration would undermine the stability of the 
current configuration. 

The CM then unambiguously identifies significant rear- 
rangements in genotype space. In what follows, monitor- 
ing the CM position allows us to explore how the tempo 
of such rearrangements evolves with the age of the sys- 
tem, counted from the beginning of the simulation.. 

Figure [T] plots the changes in this quantity, with a plot 
of the total population for reference. Importantly we see 
that both variables exhibit the same behaviour, and it 
follows that one can be used as a good proxy for the 
other. In what follows we analyse population increments 
as a proxy for the underlying shift in configuration, since 
this is much more computationally efficient. 

B. Aging 

There are two distinct tempos in the dynamics of the 
model, as has been suggested in previous papers. The 
nature of these two tempos is investigated here. Figure 
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FIG. 1: The time evolution for a typical run of the Tangled 
Nature model. The diversity (blue) displays intermittent be- 
haviour, exhibiting long plateaux interspersed with chaotic 
fluctuations (the diversity counts the number of species in 
the system, and behaves in much the same way as the pop- 
ulation). The change in centre of the mass of the system 
(red) follows the same pattern - in general there is very little 
movement, with occasional big leaps in position. 



Q] gives a qualitative picture of their origin - we see that 
the dynamics can be broadly divided into chaotic and 
stable phases, with the stable phases tending to be long- 
lived, as opposed to the chaotic phases which tend to be 
brief interludes. Within the stable phases there is still 
some dynamics due to the stochastic nature of the sys- 
tem: fluctuations and sometimes natural cycles around 
a fixed point are observed. Thus, there are two obvious 
timescales in the system - one governing the fluctuations 
within each stable phase, the other the chaotic transitions 
between stable phases. 

To determine the tempo of these phases, or more pre- 
cisely to determine how the frequency of these two types 
of events scale with the age of the system we use a tech- 
nique borrowed from analysis of spin glasses [lU, [l3j . In 
spin glasses the thermally activated evolution can be an- 
alyzed by constructing the Probability Density Function 
(PDF) of the heat flow between the system and its ther- 
mal reservoir. Large, irreversible flows of heat out of the 
system only occur when the system undergoes a signifi- 
cant rearrangement, or quake. Consequently, observing 
how the PDF evolves over time reveals how the tempo of 
reversible and irreversible events evolve over time. In 
our case, an exact equivalent of the energy in a spin 
glass is not available. Instead, we choose a character- 
istic scalar function of the dynamics, the population. In 
spin glasses the energy has a much more well defined role 
in the dynamics than does the population in the Tan- 
gled Nature model -possible changes in configuration are 
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namely accepted or rejected based on energetic consider- 
ations, whereas the population plays no such discrimina- 
tory role in the present model. 

Letting t w be the age of the system, we consider the 
distribution of population increments 

N(t w + kAt) - N{t w + (k - l)Ai), 

where k = 1, 2,3, xffe, i- e - a t different epochs in the 
system's evolution over the observation interval [t w , t w + 
%■], . The nature of the PDF will be determined by 
three parameters - t w , K and At. The first quantity, 
t w , measures as mentioned the age of the system, while 
K is inversely proportional to the size of the observation 
window. The last parameter is the number of Monte 
Carlo sweeps between consecutive measurements of the 
population size. We need the observation window to be 
large enough to capture both the quasi-equilibrium and 
non-equilibrium dynamics, while not being so long that 
the system changes too significantly over the window. In 
practice this is a matter of trial and error. In figures [2] 
and [3] the PDF of the population fluctuations is plotted 
for three different epochs, with t w = 1000, 10000, 100000 
and K = 4. The data comes from an ensemble of 1500 
runs. 

Figure [5] shows the data with At = 10 for all 
epochs. The first thing to note is that the PDF has the 
characteristic shape common to similar spin glass stud- 
ies, with one important difference. Since the data is plot- 
ted on linear-log axes, the parabolic shape in the centre 
represents a gaussian distribution, and the straight tails 
are exponential. The difference from the spin glass case 
is that here we have both positive and negative tails, 
whereas a in the heat flow PDF for a spin glass the pos- 
itive tail is undetectable, corresponding to the fact that 
quakes reduce the energy of the system. The population, 
on the other hand, may go up or down in a quake, al- 
though, as we shall see, there is still an overall tendency 
for the population to increase. 

In analogy with spin glasses, we identify the large, rare 
events in the tails with the irreversible quakes in the 
system. We see that as the system ages, these quakes 
become less and less frequent while the gaussian peak 
remains unchanged. We quantify this observation by 
rescaling the data - we anticipate a logarithmic decrease 
in the quake frequency, and so we make At proportional 
to the age of the system. When we do this, the PDF 
from three different ages collapses on to a single curve, 
so that the probability of large fluctuations (that is fluc- 
tuations appearing in the exponential tails) must scale as 
t~ , since the probability of a large fluctuation in time 
At is to leading order AtP(t w ), where P(t w ) is the prob- 
ability of a large fluctuation as a function of t w . 

The waiting time distribution for quakes (figure 0]) 
shows that the ratio of waiting time to age of the sys- 
tem remains non zero right throughout a million genera- 
tion run, giving further proof that the quake probability 
decays with the age of the system. 
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FIG. 2: The PDF of population increments for an ensemble 
of 1500 runs of 10 6 generations. Black: t m = 1000, red t w = 
10000, blue t w = 100000. In this plot At = 10. The older 
the system, the less probable a large population fluctuation 
becomes. 



IV. DISCUSSION 

We have shown that by studying the tempo of the evo- 
lution in a schematic model we can infer a great deal 
about the mode - broadly speaking we find that we can 
split the dynamics into two bits, one quasi-equilibrium 
and the other non-equilibrium relaxation. Here we would 
like to suggest a mechanism for the observed slowing 
down of the dynamics at the macroscopic scale. First, let 
us recall the work done in spin glasses on similar aging 
phenomena. It has been found there that many micro- 
scopically different systems can be modelled by a unifying 
mesoscopic principle of barrier hopping. In this picture 
the system becomes increasing stable against quakes as 
each quake corresponds to the jumping of a barrier in the 
energy landscape which is marginally larger than the bar- 
rier corresponding to the previous quake. Consequently, 
the evolution of the system is well described by record 
dynamics - each new maximum energy leads to a change 
in configuration and a new energy minimum. The energy 
landscape is such (for more detailed discussion see jl4l ]) 
that each new well is deeper than the previous one, so 
the relaxation of the system slows as records become less 
and less frequent. 

We propose that the similar relaxation dynamics seen 
in the Tangled Nature model is due to a similar scenario, 
with some caveats. One way of interpreting the spin glass 
dynamics is that the system falls into wells of increas- 
ing depth, from which it becomes increasingly hard to 
escape since the magnitude of the fluctuations of the en- 
ergy about each new fixed point remains constant. In 
the tangled nature model, it is harder to justify an en- 
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3: The PDF of population increments for an ensemble 

,6 



10000, blue tw = 100000. The total number of data points is 
equal for each age of the system; that is the sampling period 
At — scales with the age of the system. The data collapses 
onto a single curve, showing that the probability of a large 
fluctuation scales like j- 



ergy landscape picture, but we can measure the magni- 
i tude of fluctuations of the birthrate, which is a measure 
- of the systems ability to 'break out' from each equilib- 
rium state. To do this we have determined the standard 
; deviation of the birth rate as a function of the age of the 
■ system, sampled during the different equilibrium stages 
] of the dynamics. Indeed we find that there is a signifi- 
: cant decrease in the fluctuation width over a simulation 
of 10 6 generations. 

To complete the circle, we also need to understand 
why a decrease in the birthrate probability fluctuations 
1 in turn leads to greater stability. A spike in the repro- 
: duction probability of the system has an effect beyond 
simple increase in numbers - it also increases the like- 
\ lihood that new species in different parts of genotype 
space, generated by mutations, will interact and thus 
2000 destabilise the current configuration. It should be em- 
phasised that the difference between the fluctuations in 
the average probability and the fluctuations of the popu- 
lation amount to a difference in the weighting of the core 
and the periphery of the system - that is between be- 
tween vital, heavily populated genotypes and transient, 
sparsely populated ones. When considering the birthrate, 
each site is weighted equally - so the fact that we see a de- 
crease in these fluctuations while seeing a constant stan- 
dard deviation in the population is not contradictory - it 
merely shows that the fluctuations are increasingly lim- 
ited to sites, which are gradually becoming more similar 
in terms of the network of interactions Jjj surrounding 
the sites carrying the larger proportion of the population. 
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FIG. 4: Plot of a sample of quake waiting time ratios, where 
tk is the time of the k th quake. A Poisson distribution of 
quakes would give a more or less zero for long times, since 
the ratio would decay as xr for some average rate A. Instead 
we see many non zero points for long times. The clustering 
of points along the t axis is due to small 'secondary quakes' 
that are hard to distinguish from primary quakes. 



FIG. 5: The standard deviation of birth rate fluctuations de- 
creases over time. The window index, x, gives the start and 
end points of each sample window using the following conver- 
sion: t start = 10 • 2 X and t enl i = 18 ■ 2 X . Error bars show that 
while the effect is small, it is significant. 

Figure [5] shows the decrease in the standard deviation 
of these fluctuations plotted over logarithmically spaced 
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time intervals. The fluctuations were only measured in 
the metastable parts of the dynamics. This suggests that 
the system evolves in such a way that barrier hopping 
becomes less and less likely, not because of an increase in 
barrier height, but due to a decrease in the width of the 
fluctuations in the system. 

This is not the end of the story, since these fluctuations 



are themselves governed by the nature of the network 
of interactions. Why the systems chooses configurations 
less susceptible to fluctuations has not been answered in 
this paper, and will form the focus of future work. 
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